/*
Table_B01. Productivity and sentiment: Weather decomposition
Table_B02. Direct effect of rainfall on productivity
*/

* Table_B01. Productivity and sentiment: Weather decomposition
	clear all
	use "data/final_data_set.dta"
	xtset id year

	merge m:m cou year using "data/weather.dta"
	drop _merge

	duplicates drop id year, force
	xtset id year

	// generate precipitation variable that excludes extreme rain days
		gen PRCP1 = PRCP - DP1X

	// generate rates of growth and dummies
		foreach x in DP01 DP10 DP1X PRCP PRCP1 DP01_i DP10_i DP1X_i PRCP_i {
			gen g_`x' = ln(`x')-ln(l.`x')
			gen dum_`x' = (d.`x' < 0)
			replace dum_`x' = . if pd == .
		}

	// generate alternative dummy using moving average of past three years	
		foreach x in PRCP PRCP1 {
		rangestat (mean) `x', interval(year -4 -1) by(id)
		gen dum_`x'_mean = (`x' - `x'_mean < 0)
		replace dum_`x' = . if pd == .
		}
		// correlations across dummies
			pwcorr dum_P*, sig	
		
	// exogeneity assumption
		pwcorr rtfpna f_rtfpna zcci dum_PRCP dum_PRCP1 DP1X PRCP PRCP1 dum_DP1X, sig
		xtreg rtfpna zcci dum_PRCP i.year if pd != ., robust fe // choice: dum_PRCP dum_PRCP1 DP1X PRCP PRCP1 dum_DP1X 
		reg rtfpna zcci dum_PRCP
		
	// auxiliary regression
		xtreg cci dum_PRCP i.year if pd != ., robust fe // choice: dum_PRCP dum_PRCP1 DP1X PRCP PRCP1 dum_DP1X 
		predict cci_w, xb
		predict cci_nw, r
		* NOTE 1: the precipitation dummy allows for comparisons across countries as it abstracts from levels
		* NOTE 2: to use three-year moving average, replace with "dum_PRCP_mean"
		// other version: (does not allow for "predict" though)
		regxfe cci dum_PRCP if pd != ., fe(year id) cluster(id) // choice: dum_PRCP dum_PRCP1 DP1X PRCP PRCP1 dum_DP1X 

	// standardized
		foreach x in cci_w cci_nw {
		egen `x'mean = mean(`x')
		egen `x'sd = sd(`x')
		gen z`x' = (`x'-`x'mean) / `x'sd
		drop `x'sd `x'mean
		}

	xtset id year
	generate newid = id
	tsset newid year

	// weather- v. non-weather-related sentiment
		// full sample
			foreach x in f_rtfpna f2_rtfpna f3_rtfpna f4_rtfpna {
				qui: bootstrap, rep(200) seed(123) cluster(id) idcluster(newid) nodots nowarn: ///
				xtreg `x' zcci_w zcci_nw zd_* zF* if pd != ., fe
					eststo
			}
			esttab using "Table_B01_a.rtf", replace cells(b(star fmt(%9.4f)) t(fmt(%9.2f))) style(fixed) title(Table_B01) ///
			t starlevels(* .10 ** .05 *** .01) b(4) scalars(r2)     ///
			drop()
			estimates clear
		// G7 countries
			foreach x in f_rtfpna f2_rtfpna f3_rtfpna f4_rtfpna {
				qui: bootstrap, rep(200) seed(123) cluster(id) idcluster(newid) nodots nowarn: ///
				xtreg `x' zcci_w zcci_nw zd_* zF* if pd != . & dum_G7 == 1, fe
					eststo
			}
			esttab using "Table_B01_b.rtf", replace cells(b(star fmt(%9.4f)) t(fmt(%9.2f))) style(fixed) title(Table_B01) ///
			t starlevels(* .10 ** .05 *** .01) b(4) scalars(r2)     ///
			drop()
			estimates clear
		// Non-G7 countries
			foreach x in f_rtfpna f2_rtfpna f3_rtfpna f4_rtfpna {
				qui: bootstrap, rep(200) seed(123) cluster(id) idcluster(newid) nodots nowarn: ///
				xtreg `x' zcci_w zcci_nw zd_* zF* if pd != . & dum_G7 == 0, fe
					eststo
			}
			esttab using "Table_B01_c.rtf", replace cells(b(star fmt(%9.4f)) t(fmt(%9.2f))) style(fixed) title(Table_B01) ///
			t starlevels(* .10 ** .05 *** .01) b(4) scalars(r2)     ///
			drop()
			estimates clear

	// weak instrument robust tests and confidence sets for linear IV (fixed effects)
		// Cragg-Donald Wald F statistics and Kleibergen-Paap rk Wald F statistics for IV relevance, with Stock and Yogo (2002, 2005) critical values
			xtivreg2 f_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			xtivreg2 f2_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			xtivreg2 f3_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			xtivreg2 f4_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
		// AR and Wald tests
			weakiv xtivreg2 f_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			weakiv xtivreg2 f2_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			weakiv xtivreg2 f3_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
			weakiv xtivreg2 f4_rtfpna zd_* zF* (zcci = zcci_w), fe cluster(id)
		

* Table_B02. Direct effect of rainfall on productivity
	// import data
		// agriculture
			/* Legend:
				1 = thousand hectares
				2 = thousand tonnes
				3 = tonnes per hectare
			*/
			foreach a in rice maize soybean wheat {
			foreach x in 1 2 3 {
				clear all
				import excel "orig/OECD/agriculture.xlsx", sheet("`a'`x'") firstrow
				destring cou, replace
				save "data/`a'`x'.dta", replace
			}
			}
		// missed days
			clear all
			import excel "orig/OECD/missed_days_at_work.xlsx", sheet("`a'`x'") firstrow
			destring cou, replace
			save "data/missed_days_at_work.dta", replace		
	// merge with main data set
		// agriculture
			clear all
			use "data/final_data_set.dta"
			xtset id year
			foreach a in rice maize soybean wheat {
			foreach x in 1 2 3 {
				merge m:m cou year using "data/`a'`x'.dta", keep(match master)
				drop _merge
				erase "data/`a'`x'.dta"
			}
			}
			xtset id year
		// missed days
			foreach x in missed_days_at_work {
			merge m:m cou year using "data/`x'.dta", keep(match master)
			drop _merge
			erase "data/`x'.dta"
			duplicates drop id year, force
			xtset id year
		}
	// generate agricultural output measures
		// scaled
			foreach a in rice maize soybean wheat {
			foreach x in 1 2 3 {
				gen `a'`x'_rgdpna = `a'`x' / rgdpna
			}
			}
		// aggregate
			foreach x in 1 2 3 {
				gen agr`x' = rice`x' + maize`x' + soybean`x' + wheat`x'
				gen agr`x'_rgdpna = ln(agr`x' / rgdpna)
				gen agr`x'_1000 = agr`x' / 1000 
			}
	// generate scaled absenteeism measure
		gen missed_rgdpna = ln(missed / rgdpna)
	// summary statistics
		// agriculture
			// Panel A. Full sample
				tabstat agr1_rgdpna agr2_rgdpna agr3_rgdpna agr1_1000 agr2_1000 agr3_1000, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)
			// Panel B. G7 countries
				tabstat agr1_rgdpna agr2_rgdpna agr3_rgdpna agr1_1000 agr2_1000 agr3_1000 if dum_G7 == 1, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)		
			// Panel C. Non-G7 countries
				tabstat agr1_rgdpna agr2_rgdpna agr3_rgdpna agr1_1000 agr2_1000 agr3_1000 if dum_G7 == 0, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)		
		// missed days
			// Panel A. Full sample
				tabstat missed missed_rgdpna, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)
			// Panel B. G7 countries
				tabstat missed missed_rgdpna if dum_G7 == 1, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)		
			// Panel C. Non-G7 countries
				tabstat missed missed_rgdpna if dum_G7 == 0, ///
				columns(statistics) statistics(mean sd p25 p50 p75) format(%9.4f)		
	// merge with weather data
		merge m:m cou year using "data/weather.dta"
		drop _merge
		duplicates drop id year, force
		xtset id year
	// generate precipitation variable that excludes extreme rain days
		gen PRCP1 = PRCP - DP1X
	// generate rates of growth and dummies
		foreach x in DP01 DP10 DP1X PRCP PRCP1 DP01_i DP10_i DP1X_i PRCP_i {
			gen g_`x' = ln(`x')-ln(l.`x')
			gen dum_`x' = (d.`x' < 0)
			* replace dum_`x' = . if d.`x' == . // dubious command, creates issues below (missing obs)
			replace dum_`x' = . if pd == .
		}
	// generate alternative dummy using moving average of past three years	
		foreach x in PRCP PRCP1 {
		rangestat (mean) `x', interval(year -4 -1) by(id)
		gen dum_`x'_mean = (`x' - `x'_mean < 0)
		replace dum_`x' = . if pd == .
		}
		// correlations across dummies
			pwcorr dum_P*, sig	
	// dummy interaction
		gen dum_PRCP_dum_G7 = dum_PRCP * dum_G7
	// samples
		// agriculture
			tabulate country if dum_PRCP != . & agr1_rgdpna != .	
				********* Australia, Canada, New Zealand, Switzerland, United Kingdom, United States
			tabstat year if dum_PRCP != . & agr1_rgdpna != ., statistics(min max)
				********* sample starts in 1990
		// missed days
			tabulate country if dum_PRCP != . & missed_rgdpna != .
				********* all except New Zealand
			tabstat year if dum_PRCP != . & missed_rgdpna != ., statistics(min max)
				********* sample starts in 2010				
	// regressions
		// agriculture
			qui: regxfe agr1_rgdpna dum_PRCP, robust fe(year) cluster(id)
				eststo
			qui: regxfe agr1_rgdpna dum_G7 if dum_PRCP != ., robust fe(year) cluster(id)
				eststo
			qui: regxfe agr1_rgdpna dum_PRCP dum_G7, robust fe(year) cluster(id)
				eststo
			qui: regxfe agr1_rgdpna dum_PRCP dum_G7 dum_PRCP_dum_G7, robust fe(year) cluster(id)
				eststo
			esttab using "Table_B02_a.rtf", replace cells(b(star fmt(%9.4f)) t(fmt(%9.2f))) style(fixed) title(Table_TBD) ///
			t starlevels(* .10 ** .05 *** .01) b(4) scalars(r2o)     ///
			drop()
			estimates clear
		// missed days
			qui: regxfe missed_rgdpna dum_PRCP, robust fe(year) cluster(id)
				eststo
			qui: regxfe missed_rgdpna dum_G7 if dum_PRCP != ., robust fe(year) cluster(id)
				eststo
			qui: regxfe missed_rgdpna dum_PRCP dum_G7, robust fe(year) cluster(id)
				eststo
			qui: regxfe missed_rgdpna dum_PRCP dum_G7 dum_PRCP_dum_G7, robust fe(year) cluster(id)
				eststo
			esttab using "Table_B02_b.rtf", replace cells(b(star fmt(%9.4f)) t(fmt(%9.2f))) style(fixed) title(Table_TBD) ///
			t starlevels(* .10 ** .05 *** .01) b(4) scalars(r2o)     ///
			drop()
			estimates clear

			
